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Abstract 

We study the thermodynamic behaviom of an ideal gas of bosons trapped in a three- 
dimensional anisotropic harmonic oscillator potential. The condensate fraction as well 
as the specific heat is calculated using the Euler-Maclaurin approximation. For a finite 
number of particles there is no phase transition, but there is a well defined temperature 
at which the condensation starts. We also consider condensation in lower dimensions, 
and for one-dimensional systems we discuss the dependence of the condensate fraction 
and heat capacity on the ensemble used. 



1 Introduction 

The observation of Bose Einstein condensation (BEC) in ultracold alkali gases has led to in- 
tensive experimental and theoretical investigations of the phenomenon . Since the atomic 
gases of the experiments are dilute, and hence weakly interacting, their physical properties 
may be explored theoretically in terms of basic interatomic interactions and as a first ap- 
proximation the interactions may be neglected 0. This is strongly opposed to the physical 
properties of related phenomena such as superconductivity in metals and superfluidity in 
liquid He, which are fundamentally modified by strong interparticle interactions. 

In the remarkable pioneering experiment ^''Rb vapour was magnetically trapped and 
evaporatively cooled to approximately 170 nK, where a finite fraction of the particles was 
observed to occupy the ground state. When further lowering the temperature this frac- 
tion increased abruptly, signalling a Bose-Einstein condensation. The confining magnetic 
field was a time-averaged, orbiting potential |^ whose effective potential to first order is a 
three-dimensional anisotropic harmonic oscillator. The potential was not anisotropic in all 
three directions, but was cylindrically symmetric. Shortly after, a similar experiment was 
performed with ^^Na gas 0; in this case the potential was completely anisotropic. 

With these experiments in mind, we will in the present paper extend our previous work 
on the thermodynamics of an ideal isotropic gas 0, |^ to the anisotropic case. 

Although the interactions among the atoms are weak, they are quite important in the 
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condensate. Both experiments and theory show that due to the interaction the actual ground 
state is considerably larger than the single particle ground state. However, close to the 
transition temperature just a few particles occupy the ground state and the gas is much more 
dilute. At these temperatures one would believe that the ideal gas is a good approximation. 

In Ref. it has been proposed that highly anisotropic traps may freeze out one or two 
dimensions, so that the system to a good approximation is lower dimensional. This implies 
that it is of interest to study these systems also, and we present a preliminary analysis of 
them in this work. 

The plan of the article is as follows. In the next section we study noninteracting bosons 
in a confining anisotropic harmonic potential in three dimensions (3D). In the framework of 
the grand canonical ensemble we calculate the condensate fraction, the internal energy and 
the heat capacity analytically using Euler-Maclaurin summation and compare with exact 
numerical summation. In section three we discuss the same system in one and two spatial 
dimensions. In particular, we study the effect of using different ensembles in one- dimensional 
systems. In section four we summarize and conclude. 



2 Bose-Einstein Condensation in Three Dimensions 

In this section we will discuss Bose-Einstein condensation in three-dimensional harmonic 
traps. Firstly, we consider the isotropic potential and then we move on to the anisotropic 
trap. However, we would first like to discuss Bose-Einstein condensation and phase transi- 
tions in a more general setting, and comment upon earlier work. 

Let us begin our discussion by considering free bosons in 3d in the thermodynamic limit 
{N,V —>■ oo, N/V = p finite). In this case one may show (see any standard textbook in 
statistical mechanics) that the fraction of particles in the ground state is given by 

^ = l-f-)''' (1) 



where Tc is defined through 



As can be seen from the above equation, Tc defines the onset of the condensation of particles 
into the ground state. Moreover, it can be demonstrated that the chemical potential vanishes 
for T < Tc and that /i < for T > Tc. Thus, one can define the critical temperature to be 
the highest temperature for which /i = 0. The derivative of /i with respect to temperature 
is discontinuous at T = Tc, and so there is a first order phase transition. 
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Generally, the chemical potential must satisfy /i < Eq, where Eq is the ground state en- 
ergy of the system, in order to have positive occupation numbers. In p|, Kirsten and Toms 
propose that a necessary criterion for Bose-Einstein condensation is that /i must reach the 
ground state energy Eq at finite temperature. From this, they show that BEG only occurs 
for free bosons when Z) > 3,Q and not for free bosons in a finite volume in three dimensions. 
In a system with a constant magnetic field it only happens in five or more dimensions. This 
latter result has also been verified independently by Elmfors et a/ 

From theoretical work as well as from an experimental point of view we will argue that 
this criterion is somewhat too strict. In Ref. 0] it has been demonstrated numerically that 
the chemical potential never reaches its critical value in the case of bosons in a confining 
harmonic potential. According to the criterion mentioned above, this system does not ex- 
hibit BEG. However, the occupation number of the ground state changes considerably in a 
narrow temperature range, and attains a macroscopic value. More specifically, it has been 
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shown (see e.g. [|]) that the condensate near Tq = f^(^^) varies as 

^.i_(i:)%o(iv-./3). ,3) 

This abrupt change in number of particles in the ground state takes place near the temper- 
ature at which the heat capacity has a pronounced maximum (which is Tq). We would not 
hesitate to call this a condensation, although there is no critical behaviour and so it is not 
obvious how to define a critical temperature. However, in Ref. it has been proposed 
that the critical temperature is defined by the maximum of the specific heat. According to 
this criterion there is a phase transition in two and three dimensions, but not in one (see 
also section three). We think is it natural to call the temperature Tq the condensation tem- 
perature, since this temperature determines the onset of condensation of particles into the 
ground state. This view is supported by experiment, namely the fact that one has measured 
a large increase in the number of particles in the ground state over a narrow temperature 
region. 



2.1 Isotropic Potential 

The eigenfunctions of the three-dimensional harmonic oscillator Hamiltonian can always be 
written as products of three eigenfunctions of the one-dimensional oscillator, which are Her- 
mite polynomials. Hence, the eigenfunctions can be labeled by three integers n^, Uy and 
riz and the energy is Tiojijix + Uy + nz) in the isotropic case (we have set the ground state 
energy to zero). Each energy level can thus be characterized by a single quantum number 
n and the corresponding degeneracy is gn = {f^ + 1)('^ + 2)/2, which simply is the number 
of ways of writing n as the sum of the three positive integers n^, Uy and n^. The quantum 
number n then takes the values n = 0, 1, 2, . . .. Due to the high degree of symmetry of the 

^Or more precisely that Tc = when D < 2. 
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isotropic harmonic well, there is another convenient basis set. The Hamiltonian commutes 
with U' and and hence the eigenfunctions can be labeled by L, and the energy. The 
eigenf unctions are Laguerre polynomials times a spherical harmonic. 

At temperature T and chemical potential the particle number N is given by 

°° a 



n=0 



6' 



where (3 = l/ksT and gn/{e^^^" ^■^ — 1) is the usual Bose- Einstein distribution function. 
Introducing the fugacity A = exp(/3/i), the ground state particle number may be written as 

N„ = (5) 



The number of particles in excited states is then 

gnXe' 



-bn 



n=l ^ 



which is a function of the effective fugacity A = Ae^^ where h = jShuj. 



(6) 



Since the particle number is given, one must express the chemical potential as a function 
of temperature and particle number by inverting Eq. @) in order to derive thermodynamic 
quantities describing the system. Of particular interest to us is the particle number in the 
ground state as well as the specific heat of the system as a function of temperature. 

Now let us move on to actual calculations and consider different ways of handling Eq. . 
At high temperature 6 <^ 1, the difference between the terms in the sum is small, and one 
can approximate the sum by an integral. The dominating contribution at high T is found 



by setting gn — "n. /2, leading to the semiclassical limit & 11]. The transition temperature 



as defined by the temperature at which A = 1, is in this approximation, found to be 

^ huj ( N V'^ 

^0 = 7- • (7) 



In the same approximation the specific heat Cy exhibits a discontinuity at this temperature 
with the shape of a A. The discontinuity is however an artifact of the approximation. An 
exact numerical calculation 0] reveals that Cy is a continuous function of T. This is a 
consequence of the finiteness of the particle number. 

As discussed in 0, an efficient way of performing summations is given by the Euler- 
Maclaurin summation method [Q. The sum is then written as an integral plus an infinite 
series of correction terms 

n=b „{, 1 1 

E f{n) = dxfix) + -Ifib) + fia)] + -[/'(&) - f'ia)] + ■ ■ • (8) 
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Applying this formula to the sum over excited states, the contributions from the upper limit 
vanish (it is natural to single out the ground state from the sum and use Euler-Maclaurin 
on the rest. Moreover, the series converges much faster doing so). Higher order derivatives 
mainly bring down higher powers of hu/ksT, although differentiation of the degeneracy 
factor yields small corrections to lower order. At high temperature we may therefore safely 
truncate the series. Including only the first derivative, the number of excited particles 
becomes 

= Gs{X) + Ig2{X) + Gi(A) + + ^) . (9) 

We have here introduced the functions 

These functions satisfy a simple recursion relation 



G,+,(A) = -/ rfx^^^. (10) 



G,+i(A) = -— ln(l-A) + -y^ yG,(A). (11) 
They can be expressed in terms of the more familiar polylogarithmic functions 

oo n 

= E (12) 

n=l "' 

which are generalizations of the Riemann (^-function. For p = 0,1, and 2 we have 

Gi(A) = iLii(A) 
G2(A) = ^Lii(A) + pLi2(A) 

G3(A) = lLii(A) + iLi2(A) + ^Li3(A). (13) 

The Li-functions are obtained when the lower limit of the integral in Eq. ([I0|) is set to 
zero. For n < 1 they may be written in terms of more familiar functions, e.g. Lii(2;) = 
— ln(l — z), LiqIz) = z/{1 — z). The polylogarithms satisfy the functional relation 

^— = -U,^,iz). (14) 
In terms of these functions, the total number of particles reads 

Ne = ^ + ^Lii(A) + ^Li2(A) + ^Lii(A) 
1 A /31 b \ 

The internal energy as well as the heat capacity have been derived in a similar way in Ref . @ . 
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2.2 Anisotropic Potential 

The potentials used in the experiments 0, Q are either partially or completely asymmetric. 
We therefore need expressions for the cases where the frequencies are different. Let us first 
consider the case uj^ = '^y '^z- As for the isotropic potential we use the Euler-Maclaurin 
formula for each of the sums. The calculations are similar and we skip the details. The 
result is 



1 2 11 

+— Li2(A^.) + — Lii(A^.) + — Lio(A^.) 

+^Lii(A,) + ^Lio(A,) + 0(Li_i). (16) 

We have here introduced the notation X^z = Aexp(— 6^. — b^) etc. where = hjSuJx and so 
on. We also need the internal energy and the result using the Euler-Maclaurin formula is 

U,T = -i-L,(A„)+(^ + l)u3(A„) 

(7b, 7 35 \ [ bl 5b, 23 176^.\ 

+ 1,1261 + 26: + 126: j + 1,1261 + T26: + T2 + T26:j 

( bl 796, 196, bl bl \ 

+ 1,66: + 144+^^ 36061 +126: 

2 3 56 

+ ^Li3(A,) + ^Li2(A,) + 2Lii(A,) + ^Lio(A,) 

1 ^6 
+ -Li2(A,) + Lii(A,) + -fLio(A,) + 0(Li_i). (17) 
b,. 12 



The specific heat is now found by differentiating U with respect to T. In doing so we 
need dfi/dT — fJ^/T which is found by differentiation of A^. 

In the ^^Na experiment all three frequencies were different [0 . Thus, we need the expression 
for the particle number in this case also. By applying the Euler-Maclaurin approximation 
one finds that reads 
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Here, "sym." means symmetrization, e.g. 
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The two leading terms at high temperature are 
N 



1 T- /XN 1 

Li3(A) + - I ^ + 
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(20) 



In [jT^ Grossmann and Holthaus considered a potential where ujy = \/2ujx and uj-, = \/^ijJx- 
They constructed a continuous density of states based upon the two leading terms in the 
degeneracy. For the Li2-term they found a coefficient 'j{kBT/hu) where u = (cu^^cUyCU^)^^'^- 
A numerical summation gave 7 ~ 1.6. This numerical factor is now easily found from Eq. 

my- 



7 



3\V3 / 1 1 1 \ 



(21) 



which agrees quite well with their approximate result. The same value for 7 has been ob- 
tained independently in |14]. 



For the internal energy one finds 
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+Lii(A^) +sym. + C(Lio). 
Again the specific heat can be found by differentiation of this expression. 



(22) 
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2.3 Results 



With the formulas derived above at hand, we can make predictions for the actual physical 
systems within the ideal gas approximation. In doing so we have used the values for the 
oscillator frequencies and particle numbers as given in [0, ^. 

For the completely anisotropic system (^^Na) the frequencies are 410, 235 and 745 Hz, and 
the particle number N is 700 000. The frequencies for the ^''Rb experiment are Uz = 2407r 
Hz and u!x,y = ^z/V^, and the particle number is 20 000. For these systems we have 
calculated both the condensate fraction and the specific heat using the Euler-Maclaurin 
approximation. The results are compared with those of direct numerical summation and 
the semiclassical approximation. The temperature is measured in units of the semiclassical 
transition temperature Tq. It is 1666 nK for the sodium gas and 127 nK for the rubidium gas. 

In Fig. 1 we plot the condensate fraction as function of temperature. The numerical and 
analytical results agree well. For N = 700 000 the transition is quite sharp, and the conden- 
sate nearly vanishes about 1% below Tq. For the smaller particle number, the transition is 
smoother. 



0.50 




Figure 1: The condensate fraction Nq/N in three dimensions with anisotropy-parameters taken from the 
^^Na and ^^Rb experiments. Nq/N is plotted as a function of T /Tq. 

The specific heat is shown in Fig. 2. We again find good agreement between numerical 
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and analytical results. The discontinuity in the semiclassical approximation is smoothened 
by the corrections. Still, for = 700 000, the specific heat falls rapidly over a temperature 
range of less than 10 nK. For this case it has a maximum value of Cy/Nks — 10.43 , quite 
close to the semiclassical value 12^(4)/(^(3) ~ 10.80. 
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Figure 2: The specific heat Cy in an ansiotropic three-dimensional trap with frequencies taken fr'om the 
^^Na experiment. Cy is plotted as a fmiction of T/nK. 



3 BEC in Lower Dimensions 



All BEC-experiments on trapped Bose gases reported so far are performed on three-dimensional 
systems. However, by the use of optical dipole traps or highly anisotropic magnetic traps 
one might, effectively, freeze out oscillations in one or two of the dimensions at low temper- 
atures. In Ref. [0 the authors suggest the use of a field configuration with a radial energy 
level spacing of 200 nK, and a much smaller axial spacing. At temperatures below 100 nK 
this system should effectively be one-dimensional. 

In the previous section we discussed in some detail different criteria for BEC. The same 
arguments apply here and we make a few comments. For a free gas in the thermodynamic 
limit, the chemical potential vanishes first at T = 0. This implies that the critical tempera- 
ture in one and two dimensions is zero and that the condensate vanishes at finite temperature. 
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In a finite volume there is again no critical behaviour and there is a condensate even at finite 
temperature. This is also the case for an ideal gas trapped in an external harmonic potential, 
as we shall see below. 

In this section we consider the trapped Bose gas in both one and two dimensions using the 
same methods as for three dimensions. We start by examining the case of a highly anisotropic 
three-dimensional system, and demonstrate the effective reduction to one dimension when 
two of the oscillator frequencies are much higher than the third. After this, we briefly discuss 
the use of different ensembles in one dimension. As opposed to two and three dimensions it 
is here quite straightforward to use the canonical ensemble. For completeness, we also give 
some results for the two-dimensional case. 



3.1 Reduction to One Dimension 



Excitations in a given direction of the trap will be suppressed when the temperature is below 
the energy scale set by the corresponding oscillator strength. The highly anisotropic trap 
with ujx,y ^ should therefore effectively be one-dimensional at temperatures ksT <C 
hujx,y In Fig. 3 we have plotted the condensate fraction and speciflc heat Cy for the case 
^x,y = SOOuJz- In the same flgure we also give the results for a one-dimensional system with 
frequency Uz- 



The three-dimensional condensate fraction is seen to agree very well with the one- 
dimensional one at low T, but at higher temperatures the former is somewhat reduced 
due to radial excitations. These excitations have a more dramatic effect on the speciflc heat. 
The two curves have high temperature limits of three and one, respectively. 



The expressions for the particle number and the energy density in the one-dimensional 
case are given by: 



n=0 



(23) 



U 



oo 

E 

n=0 



j/3(e„-At) _ ]^ 



^U,{\) + iLii(A) + ^Lio(A) + ^Li_i(A) 



(24) 



This one-dimensional system has earlier been studied by Ketterle and van Druten 0. 
For (3uj <^ 1, they obtained the equation 

N = No + ^Lio(Ae-^'^/2). (25) 
b 

The contribution from the excited states is similar to the flrst term in Eq. ([23|), but has a 
different effective fugacity. From Eq. (p5|) , Ketterle and van Druten have deflned a transition 
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Figure 3: Comparison of the condensate fraction and the specific heat of the highly anisotropic trap and 
the one-dimensional trap. 



temperature Tq by setting A^o = and A = 1. We then have ^ -^f^ In ^^7^. In Fig. 4 we 
compare various approximations for the density, plotting the condensate fraction as function 
of temperature for N = 2 ■ 10^. Including all three terms in Eq. (23), the analytical and 
numerical results agree to three decimal places. The heat capacity is a smooth function of 
temperature in the limit — ^ 00, and so there is no phase transition in one dimension. 



3.2 Comparison of Ensembles in One Dimension 

The trapped Bose gas has so far been described in terms of the grand canonical ensemble. 
This implies an open system exchanging energy and particles with a reservoir, where the 
average particle number is given by the chemical potential. 

The experiments are performed at fixed energy and particle number, and should therefore 
be described in terms of the microcanonical ensemble. The derivation of the corresponding 
partition function is in general a difficult task, which involves complicated combinatorial 
problems. Using the canonical ensemble simplifies the calculations considerably, and still 
satisfies the constraint of fixed particle number. 



The different ensembles are known to give the same predictions for thermodynamic quan- 
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Fi gure 4: The condensate fraction in different approximations as a function oiT/uj. N = 2 ■ 10^ 



titles in the large limit |T^. Below we compare the results for the grand canonical and 
canonical ensembles. The results for the specific heat turn out to agree well. However, for 
the condensate fraction we find significant deviations. 



Assuring that each configuration is counted only once, the partition function in the 
canonical ensemble can be written as 

oo oo W -1 

^ = E ■ ■ ■ E = n T— 7- (26) 

ni=0 n]v=n]v_i i=l 

Here, x = exp{—(3uj). The internal energy is now easily found to be 

The derivation of the condensate fraction is a little more involved: It can be written as 

N 

iVo = E ^Pfc' (29) 

k=l 



and the specific heat becomes 
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where pk is the probabihty of finding k particles in the lowest energy level. This probability 
is found summing all contributing configurations: 



oo oo 

Pk = j ■■■ ^ 



nfc+i+"fc+2+---+"iV _ x'^~'' 



k-1 



na 

i=0 



(30) 



In Fig. 5 we have plotted the condensate fraction and the specific heat for N = 10^ and 
= 10^, as given by the two ensembles. The difference in specific heat is hardly visible for 
the lower particle number, and for = 10^ the results of the two different ensembles agree 
to six decimal places. 




T/T. 



Figure 5: The condensate fraction and the specific heat in different ensembles as a function of T/Tq and 
for different values of N. 



The situation for the condensate fraction is very different. One can show that 
the error made calculating A^q in the grand canonical ensemble is, at most, of the order 
0{1/ InN). This is seen from the numerical results given in Fig. 5, where a significant 
deviation is found at intermediate temperatures even for N = 10^. 

This result is in contrast to the case of a free gas in three dimensions, where the deviation 
is known to go like (lnA^)/A^ [jl8l. A more complete discussion of the finite A^ effects for this 
system in the canonical ensemble will be given elsewhere [O . 
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3.3 Two Dimensions 



For completeness we also include results for an isotropic two-dimensional trapped gas. The 
equation for the particle number here takes the form 



N = No + ^Li2(A) + ^Lii(A) + ^Lio(A) + ^Li_i(A). 
0'^ 12 b 



The expression for the internal energy reads 
U = 00 



|Li3(A) + ^Li2(A) + ^Lii(A) + hio{\) + ^Li_i(A) 



(31) 



(32) 



For large one may define an approximate condensation temperature setting A^o = and 
A = 1. The leading term at high temperature then gives [0 



Tn\2 ... Tn IN 1 



iV=(— L^2(l) or ^=/__ = ^V67V. (33) 



UJ J UJ \l C(2) TT 

The limiting value at large N for the condensate fraction is 



At temperatures just below Tq and large A^ the energy takes the value U = 2uj({3){T/ujY. 
This gives a specific heat 

C^(T^To-) = 6A^C(3)/C(2). (35) 
Above Tq the leading terms for Cy are, again in the large A^ limit 

Cv MX) MX) 

W = 'iMX) - ^iMxy ^^'^ 

Setting A = 1, we see that the second term vanishes, and the specific heat approaches the 
same value as found above. There is thus no jump in the specific heat even in the large A^ 
limit in two dimensions. 

In Fig. 6 we have plotted the condensate fraction for different values of A^. The transition 
sharpens as A^ grows. The specific heat displayed in Fig. 7 is seen to develop a sharp, though 
continuous peak with increasing A^. 



This system behaves similarly to the one in three dimensions: The specific heat has a 
maximum, and the condensate has a powerlike fall off in the large A^ limit (A"o ~ 1 — (T/Tq)^ 
in 2D and A'o ~ 1 — (T/Tq)^ in 3D). Thus, as opposed to the free case, a two-dimensional 
trapped Bose gas undergoes a phase transition to a condensed state at finite temperature^. 

^During the preparation of this paper we came across a recent preprint by W. J. MuUin [jisl , where the 
same conclusions concerning BEC in lower dimensions are drawn. 
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Figure 6: The condensate fraction for different values of the particle number TV as a function of T jT^. 



4 Conclusions 



In the present paper we have discussed the thermodynamics of an ideal Bose gas trapped in 
an anisotropic harmonic well in some detail. We have investigated the condensate fraction 
and the heat capacity for realistic values of the frequencies by numerical as well as analytical 
methods. There is an excellent agreement between the two methods in all dimensions. 
Moreover, we have studied the highly anisotropic trap where one of the frequencies is much 
smaller than the others. This system should be realizable in future experiments and it is 
expected to be effectively one-dimensional for low temperatures. We have therefore compared 
the condensate and the heat capacity for this system with the corresponding quantities in 
the one-dimensional trap. The overall agreement for low temperatures confirms the above 
expectations. We have also considered thermodynamic quantities in different ensembles at 
finite N . The most remarkable result is the deviation of the number of particles in the 
ground state. 

Finally, we would like to mention several possible improvements of the present treatment 
of BEC. The most obvious is the inclusion of interactions between the atoms in the gas. This 
could either be based on perturbation theory or some nonperturbative method. There are 
already a few articles on the subject |T^, but more work needs to be done. Secondly, one 
should use the microcanocial ensemble in order to compare with experiment. It is essentially 
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Figure 7: The specific heat for different values of the particle number as a function of T/Tq. 



a combinatorial challenge to do so. More insight into finite effects of the ideal system is 
also of interest, and a detailed discussion will be presented in a follow up paper by one of 



the present authors 
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